x <- seq(-8, 8, by=0.5)
v <- foreach(y=x, .combine="cbind") %dopar% {
  r <- sqrt(x^2 + y^2) + .Machine$double.eps
  sin(r) / r
}

persp(x, x, v)
